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Two-dimensional thermosolutal convection with no-slip boundary conditions is studied using numerical 
simulations in a periodic domain. The domain is large enough to follow the evolution of phase instabilities of 
fully nonlinear traveling waves. In the parameter regime studied these instabilities evolve, without loss of 
phase or hysteresis, into a series of confined states or pulses characterized by locally enhanced heat and solute 
transport. The wavelength and phase velocity of the traveling rolls within a pulse differ substantially from 
those in the background. The pulses drift in the same direction as the convection rolls on which they ride but 
more slowly, and are characterized by an exponential leading front and an oscillatory trailing end. Multiple, 
apparently stable, states are found for identical parameter values. The qualitative properties of the pulses are in 
good agreement with the predictions ol a third-order phase equation which accounts for the relation between 
wave number and phase velocity, the oscillatory tails and the multiplicity of states. These properties of the 
pulses are shown to be a consequence of Shil’nikov dynamics in the spatial domain. [S 1063-65 lX(98)04901-0] 
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I. THE PHYSICAL PROBLEM 

The dynamics of overstable convection produced by the 
competition between stabilizing and destabilizing mecha- 
nisms in doubly diffusive systems has attracted much inter- 
est. Typical of such systems is thermosolutal convection in 
which this instability usually develops into nonlinear travel- 
ing waves. Such waves were studied extensively [1,2] using 
a combination of bifurcation analysis and numerical simula- 
tion in two dimensions, employing periodic boundary condi- 
tions with period 2n/k c in the horizontal. Here k c is the 
wave number with which the instability sets in. While such 
simulations are suitable for the study of the relative stability 
between traveling and standing waves they exclude long- 
wavelength modes to which such waves can become un- 
stable. 

This paper is devoted to the study of long-wavelength 
instabilities of fully nonlinear traveling waves. For this pur- 
pose we employ a relatively large aspect ratio T = 64 with 
periodic boundary conditions applied at .r = 0, 64 allowing 
something of the order of 26-32 roll pairs to form in the 
spatially uniform state, depending on parameters. We find, 
and describe in detail, a new dynamical traveling wave state, 
consisting of one or more nonlinear pulses or confined states, 
comprising rolls of different wave number, phase velocity, 
and amplitude than the background wave train, that travel 
more slowly than the background while transporting more of 
the heat and solute. We provide convincing evidence that 
these pulses form as a result of a phase instability of the 
background wave train, and describe a qualitative theory 
based on a third-order phase equation that accounts for their 
structure. In particular, we suggest that both single-pulse and 
multipulse states are manifestations of Shil’nikov dynamics 
in the spatial variable and thereby account for the observa- 
tion, at fixed parameter values, of multiple, numerically 
stable, pulse states. 

Confined states have been the subject of much interest in 
the last several years, motivated largely by the discovery of 


such states in binary fluid convection [3]. This system is 
closely related to thermosolutal convection but differs from 
it in several crucial ways. In binary fluid convection the sta- 
bilizing solute gradient develops in response to the applied 
adverse temperature gradient; no solute is therefore trans- 
ported through the system. The traveling waves that are ob- 
served appear via a subcritical Hopf bifurcation, and very 
close to onset are seen to evolve into patches of traveling 
waves separated by an essentially quiescent fluid. The enve- 
lope of these waves moves rather slowly and in some experi- 
ments may be stationary. A qualitative theory of this phe- 
nomenon, based on envelope equations coupled to the 
concentration field, appears to account for most of the ob- 
served properties of these confined states. It should be noted 
that the subcriticality of the bifurcation is an essential part of 
this theory [4], which describes the traveling wave patches in 
terms of solitary waves coupled to the large scale concentra- 
tion field [5], Such waves are in turn to be viewed as ho- 
moclinic orbits to the origin in an appropriate amplitude 
equation, that is, the amplitude increases from zero at 
a = -oo to a maximum and then decreases to zero again as 
x— >+oo [6], Such states have also been found in numerical 
simulations of the full partial differential equations with ex- 
perimental parameters and are in excellent agreement with 
their measured properties [7], 

The confined states that form the subject of the present 
paper are quite different. In particular, they connect a uni- 
form but fully nonlinear wave train at x=-x to itself as 
x ~ >0 °- Consequently, the traditional amplitude equation ap- 
proach is inadmissible; the pulses must be described using a 
real order parameter, the phase of the wave. We are able to 
locate numerically a regime in which the phase of the wave 
is apparently conserved, and develop our theory based on 
phase conservation. As discussed further below, the phase is 
not always conserved and when it is not our theory no longer 
applies. 

The recent experimental realization of doubly diffusive 
convection by Predtechensky et al. [8] raises the possibility 
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that the pulses described below could be observed in an ap- 
propriately focused experiment. In these experiments the de- 
stabilizing buoyancy force is produced by an imposed con- 
centration gradient of a second species, instead of an adverse 
temperature gradient. Consequently the Lewis number (i.e., 
the ratio of the diffusivities of the stabilizing and destabiliz- 
ing components) lies typically in the range 0.1 — 1.0, a param- 
eter range that forms the focus of the present study. How- 
ever, in the following we retain the traditional thermosolutal 
terminology. 

This paper is organized as follows. In the next section we 
describe the system studied and some of its elementary prop- 
erties. In Sec. Ill, which forms the bulk of the paper, we 
describe the results of our numerical simulations. In particu- 
lar, we discuss the properties of the various fields (tempera- 
ture, concentration, vorticity) inside the pulses and outside. 
We provide detailed information on the phase velocity, wave 
number, and amplitude of the rolls inside and outside a pulse 
at three different Rayleigh numbers, as well as information 
about the speed of propagation of the pulses. We also de- 
scribe the transport properties of the pulses in both vertical 
and horizontal directions, and compare their transport effi- 
ciency with the background wave train. This information 
forms the basis of the theoretical interpretation of the results 
advanced in Sec. IV. The paper ends with a brief discussion 
of the results and an intimation of future work. 

II. THE MATHEMATICAL MODEL 

The nondimensionalized equations describing two- 
dimensional thermosolutal convection can be written in the 
form [1] 

-[^v 2 'p+y(^,v 2 'p)]=/?r^o-^' ? A+v 4 ^, (i) 

C 7 


(a) (b) 



FIG. I. Neutral stability curves (a) Rj vs k and (b) o> vs k for 
r = 64.0, R s = 30 000.0, cr = 1 .0, and r= 10" m . Dots correspond to 
the discrete modes of the finite domain, which increase in number 
in proportion to T. The continuous-case critical values provided by 
linear stability theory are /? ( r 0) = 25 251.76, £ ,0) = 4.37, o/ 0) 
= 70.66. 

The velocity boundary conditions translate into the following 
boundary conditions on the stream function: 

¥(*,0,0 = 0, ¥(*,l,0=-w(0. (6) 

<? z ¥( x , 0,0 = 0, ^¥(*,1,0 = 0, 

with m{t) determined as part of the solution. The absence of 
Galilean invariance implies that m(t) is unique; in fact m(t) 
is the horizontal mass flux resulting from any mean flow 
(£/(z,r),0) that accompanies the solution: 

r i l fi r r 

”<«- J„ Wz -»*‘-f J/'J/'ir <7) 


, 9 , 0 +•/(¥. e ) =<?,¥+ v 2 e , ( 2 ) 

<9,2 +7(^,2) (3) 

where 'R(jc,z,f) is the stream function, ®(x,z,t) and 
l(x,z,t) denote departures of the temperature T and concen- 
tration S from their linear conduction profiles in the absence 
of convection, and J(W ,/ ) denotes the Jacobian 
ft \J / j )/ft x ,z). In terms of the stream function the velocity 
is given by The dimensionless pa- 

rameters R t and R s are the thermal and solutal Rayleigh 
numbers, respectively; we shall use the former as our bifur- 
cation parameter. The quantities <j and t denote the Prandtl 
and Lewis numbers. 

Equations (1)— (3) are to be solved in the domain (x,z) 
e[0,r)X[0,l], for f>0, subject to periodic boundary con- 
ditions in the horizontal. 

(^,0.2K-r,z,f) = (^.0.S)U + r,j,t), (4) 

and no-slip, fixed temperature and concentration boundary 
conditions at top and bottom, 

„ = M , = 0 = 2 = 0, z = 0,l. (5) 


For traveling waves m(t)^ 0. 

As in earlier studies [1] we use the parameter values 
( 7 = 1 . 0 , 7-=10 _1/2 , with R s = 3X10 4 and T = 64.0. Figure 1 
shows the neutral stability curve Rf as a function of the 
horizontal wave number k and the corresponding frequency 
w(k), both for 1=64.0 and r = ®°. The dots indicate the 
discretized values corresponding to the finite aspect ratio. 
For these parameter values, the critical (thermal) Rayleigh 
number, wave number, and oscillation frequency are found 
to be 

R ( r 0) = 25 251.76, 

Jt (0) = 4.37, (8) 

tu <0) = 70.66. 

The large value of o> <0) indicates that the system is far from 
the codimension-two point. 

With periodic boundary conditions the system (l)-(5) is 
equivariant under the group 0(2) of rotations and reflections 
of a circle. As a consequence of this symmetry the number of 
eigenvalues on the imaginary axis is doubled, and two 
branches of solutions, standing waves (SW) and traveling 
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waves (TW), bifurcate simultaneously from the conduction 
solution. Specifically, near onset, 

+ 0 ( 2 ), ( 9 ) 

where A L (t) and A*(f) are the (real) amplitudes of left- and 
right-traveling waves, satisfying the equations 

A R = (\ + aA 2 L + bA 2 , 

A L =(\ + aA 2 R +bA 2 +-")A L . (10) 

Here A 2 =A 2 r + A 2 l is the square of the total amplitude, and 
\±iil are the (double) eigenvalues of the linear stability 
problem. Generically, these equations possess only three per- 
sistent solutions: the conduction state (A L , A /? ) = (0,0), left- 
traveling waves (A,0), and standing waves ( A/V2 , AiVT)\ all 
other solutions can be obtained using translations and/or re- 
flections of these elementary solutions. The stability of these 
solutions is determined by the real coefficients a , b as de- 
scribed by Knobloch [9]. With stress-free boundary condi- 
tions at top and bottom 6=0 [1,9]; this is not so, however, 
for the more realistic no-slip boundary conditions (5). Re- 
nardy [10] shows that for the parameters employed both TW 
and SW bifurcate subcritically at onset. Thus neither branch 
is stable near onset. The large-amplitude traveling waves we 
find numerically at supercritical Rayleigh numbers probably 
acquire stability at finite amplitude via a saddle-node bifur- 
cation followed by the shedding of a branch of modulated 
traveling waves, as described in Ref. [11]. 

Equations (l)-(3) with boundary conditions (4)-(6) were 
solved using a hybrid finite-difference and spectral scheme. 
Spectral decomposition in the horizontal with 16 grid points 
per unit length provided 512 Fourier modes for the spectral 
decomposition at each z level of the grid. A 33 point equal- 
length finite-difference grid in the vertical provided an eco- 
nomic compromise between the resolution needs of the 
boundary layers and those of the long-range fluxes. We 
worked with the equations of motion in their biharmonic 
form (l)-(3) rather than the more traditional stream- 
function-vorticity formulation in order to avoid (a) inversion 
of the Poisson operator, and (b) difficulties arising from the 
coupling between the stream function and vorticity [12]. 
These are associated with the peculiar nature of the boundary 
conditions for the stream function ¥ and the initial condi- 
tions on the vorticity a>, which lacks boundary conditions 
[13]. In the biharmonic formulation the specification of both 
Dirichlet and Neumann boundary conditions for ¥ poses no 
difficulty because both are required to supplement the fourth- 
order elliptic operator V 4 present in Eq. (1). 

Solutions were first computed with fictitious boundary 
conditions for the stream function ¥, namely, ¥ = <9 Z ¥ = 0, 
at z = 0,l. This allows for larger time steps in the numerical 
code and a fast approach to a stationary solution. These 
boundary conditions cancel the overall lateral fluid transport, 
but the dynamics that result from them are close to the true 
ones when the net horizontal transports are small. One simu- 
lation was run from random initial conditions until the solu- 
tion matured, and then solutions for other values of R r were 
found by perturbing the previous ones. Once stationarity was 


achieved in this manner, the resulting fields were used as 
starting points for subsequent simulations in which the real 
boundary conditions for the stream function were used 
[^\p = 0, at z = 0,1, with ¥(*,0,f) = 0, ¥(*,1,0 = -m]. 
Other solutions were run fully from random initial data for 
checking purposes. 

III. RESULTS 

The nonlinear simulations reported below' explore the 
Rayleigh number range from Rj= 1.15/?^ to 1.55/? { / ) and 
cover the regime from just prior to pulse formation to the 
regime where phase is no longer conserved and solutions 
appear weakly chaotic. We concentrate on the range charac- 
terized by phase conservation, and present detailed results 
for three values of R T , viz, 1.19/? ( /\ 1.23 and 
1 21R ( T o) . 

The efficiency of convection is measured by the thermal 
and solutal Nusselt numbers. We define these dimensionless 
numbers as the ratios of the time-averaged vertical heat or 
solute flux to the corresponding conductive flux in the ab- 
sence of convection through a given plane z = const. Specifi- 
cally, we define the thermal and solutal Nusselt number den- 
sities N t (: c,z) and N s (x,z) in dimensionless variables to be 

I dT{x,z,t) . \ /f 

N t {x,z) = ( — + w(*,z,r)7X*,z,f ) ) (11) 

and 

I dS(x,Zj) , \ 

N s (x,z) = ( — + T l w(x,zd)S(x,z,t)U 

( 12 ) 

where ( ) indicates the time average. The thermal Nusselt 
number N T (x) is then given by 

JV 7 (*) = {'N-Axrfdz. (13) 

Jo 

In the simulations reported here no qualitative difference was 
found between the values of the vertically averaged Nusselt 
number N T (x) and those of N T (x,z) evaluated at midlayer. 
The total thermal Nusselt number is usually defined as 

1 f r 

Nj=j; N T {x)dx, (14) 

1 Jo 

We shall find it convenient to also introduce reduced thermal 
Nusselt numbers, defined by 

05) 

These measure the transport properties of heat in subdomains 
A of the overall domain of extent T. Solutal Nusselt numbers 
N s (x), N s , and N St ± are similarly defined. 

For a steadily traveling pattern with or without a pulse 
(but no phase slips) all global quantities are time indepen- 
dent. This is not so, however, for the reduced Nusselt num- 
bers which are modulated because of the passage of rolls into 
and out of the pulse. In addition to the above quantities the 
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local kinetic energy density E(x,zJ) = 511-11, and the enstro- 
phy density Z(x,zJ) = \a>\ 2 are used as diagnostics, with E 
and Z being the values of the associated integrals over the 
full domain. Here g>=VXu is the vorticity. 

Lateral and vertical fluxes of heat and solute have been 
computed, as have the mean lateral velocity and mass flux. 
When the fluid motion attains a statistically stationary re- 
gime these quantities become time independent. The wave- 
length of the convective rolls both inside and outside the 
pulses is an important diagnostic of the phase and has been 
carefully measured. The propagation speeds of the rolls in- 
side and outside the pulses were computed by tracking zeros 
of the stream function 'P at midlayer. The presence of a 
“mirror-glide” symmetry [2,14-16] makes the midline the 
natural place for these measurements. Finally, moving aver- 
ages of several of the fields have been obtained. Such aver- 
ages are computed in reference frames moving at the pulse 
velocity. This procedure filters out small-scale structures, re- 
vealing the pattern of large-scale behavior. 

A. Structure of the fields 

At the lower end of our Rayleigh number regime random 
perturbations of the conduction state 'R = 0 = X = 0 evolve 
into a spatially uniform traveling wave train; left- and right- 
traveling wave trains are equally likely. Although of modest 
amplitude these solutions reveal that both 0 and X have a 
trapezoidal, highly nonharmonic shape that cannot be de- 
scribed quantitatively by a weakly nonlinear expansion near 
onset [15,16]. In these states all the fields have the “mirror- 
glide” symmetry 

f(x,zJ)= -f |jc+ -,l-z,rj, ( 16 ) 

which is a combination of a translation by half a wavelength 
in the x direction with a mirror reflection through the 
plane [2,14-16]. This symmetry persists even at large Ray- 
leigh numbers as long as the number of roll pairs is con- 
served during their time evolution. Associated with these 
states small-amplitude long-wavelength modulations were 
sometimes observed; these are believed to be long-lived tran- 
sients. 

As the Rayleigh number R T is increased propagating 
pulses on a background of traveling rolls appear. When the 
simulations start from a random perturbation of the conduc- 
tion state at these values of R T there appear several such 
domains with different wavelengths and with varying widths 
and propagation velocities. During an initial transient regime 
these domains move and recombine until they attain their 
ultimate number. When two such pulses remain after this 
first stage, their evolution continues until they reach a final 
equilibrium separation that is preserved (with slight modula- 
tions in time) as they move, both at the same speed, in the 
direction of the background rolls. Simulations carried out 
with Rayleigh numbers at 1.31 /? ( r 0) or above produced states 
that failed to settle down into a stable coexistence of the 
different domains. 

We now define the quantities that will be used to report 
our findings and describe how they were measured. For this 
purpose we split the domain into the pulse and the back- 


ground, indicated by subscripts p and b, respectively, such 
that A„+A* = r The boundary separating these regions is 
defined in terms of the wavelength as measured from 
at z=j; the pulse region is defined as the region 
with less than average wavelength, with the background re- 
gion forming the complement. 

Roll wavelength \ r (x): The length of each roll was mea- 
sured at midlayer and used to define the local semiwave- 
length. A discrete function was thus built whose ordinates 
are the midpoint coordinates of each roll and whose abscis- 
sas are twice the measured length of the corresponding roll. 
The roll wavelength \ r (x) was obtained by fitting a continu- 
ous curve through these points. The average roll wavelength 
\ r was computed from the continuous fit using the expres- 
sion 

“ 1 f r 

K=y J o ^r(x)dx. 

The roll wavelength in a pulse was computed using 

\ r p « |* X r (.r)<7x, 

S p Js p 

where the subportion S p of the pulse region involves rolls 
with nearly uniform wavelengths, namely, with wavelengths 
within 20% of the relative minimum value in the pulse (cf. 
Fig. 10). The wavelength in the background was obtained 
from 

Kt = ir f 
J\ b 

where A b is the background region of the traveling pattern. 

Roll phase velocity v r (x): A discretized function was 
computed by measuring the phase velocity of the points at 
midlayer where the stream function vanishes. A continuous 
interpolation of the sampled values yields i; r (jt). The aver- 
age roll velocity v ~ r , and the roll velocity within the pulse 
v r p and in the background v rb were then computed in an 
analogous way to their wavelength counterparts. These quan- 
tities are distinct from the pulse velocity v p . 

1. R t =1J9R ( ^ 

Figure 2(a) shows a space-time plot of a one-pulse state at 
R t = 1.19fl ( r o) . Only clockwise streamlines are shown. The 
pulse travels in the same direction as the background travel- 
ing rolls, but at a smaller speed: the rolls travel with an 
average speed of v r ^ 14.21, while the pulse moves at u p 
^5.88. The complete solution consists of 26 pairs of con- 
vening rolls with an average wavelength \ r *^2.12, but the 
local wavelength inside the pulse (\ rp ^1.31) is seen to be 
approximately one-half the wavelength in the background 
2.86). These wavelengths are almost constant, except 
in the vicinity of the interface between the pulse and the 
background. The presence of the pulse thus indicates stable 
coexistence of regions with rolls of different amplitudes, fre- 
quencies, and wave numbers. The pulse is a moving region 
that the rolls have to traverse. In order to do so they shrink 
laterally, reduce their overall horizontal speed, and spin more 
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FIG. 2. Traveling waves and pulses in a sequence of equally spaced snapshots of the stream function in the full computational 

domain, with time increasing upward. Rolls and pulses both propagate rightward, (a) One-pulse state for R r = 1.19 R { } 0) involving 26 roll 
pairs, \,^2.12, \ r ^2.86, 1.31, v r ^ 14.21, v r b ^ 15.01, v rf) ^ 10.05, v p ^5.$%. (b) Two-pulse state for R r — 1.23 Rj >] involving 26 

roll pairs, \ r ^2.38, \ r ^3.68, \ r p ^\.2\, 13.11, iT rJ) ^ 15.05, u" r ^*=8.57, ^^5.42. (c) Tw f o-pulse state for R r = 1.27 Rj }) involving 32 

roll pairs, \,~2.45, \,^3.43, X^^l.ll, iT r ~ 10.66, v r h ^ 12.83, v r ^1.24 y v p ~4A9. 


vigorously until they pass out of the pulse region. For com- 
parison, at R t = 1.15 Rj >) (not shown) the solution also con- 
sists of 26 roll pairs but the pulse is absent and the wave- 
length uniform. Careful numerical investigation indicates 
that the pulse appears once R T exceeds a well-defined thresh- 
old and that the formation process is nonhysteretic. 


Figure 3(a) shows a closeup view of the fields. The tem- 
perature contours are smoother than the solute ones; due to 
its small diffusivity, solute is transported predominantly by 
fluid advection with very little diffusion, while the high dif- 
fusivity of the temperature smoothes out its contours more 
effectively. Plume structures can be seen to form inside the 
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(o) R t = 1.19 f 4 o) 


4 x 0 36 x 40 



0 x ► 64 


(i b ) R t = 1.23 R[ o) 


14 x 18 34 x 38 



0 x ► 64 


(c) R t = 1.27 Rj o) 


8 x 12 22 x 26 



0 x ► 64 


FIG. 3. Closeup views of the field structures inside the pulse 
region (left panels) and in the background (right panels) at a par- 
ticular instant in time, showing (from bottom to top) the stream 
function A r , the temperature perturbation 0, the concentration per- 
turbation X, the temperature T , and the concentration S. (a) R T 
= \A9R l r o) . (b) R t =\. 23R { t o) . (c) R t = 1.27 /?</> . 

pulse with a characteristic mushroom shape caused by the 
lateral bending of the upward plumes of heat and solute. 
Although the stream function in Fig. 3(a) has nearly the same 
number of contour lines inside and outside the pulse the 


contour lines inside the pulse are more tightly packed, indi- 
cating increased gradients, particularly in the horizontal, and 
hence increased vertical velocities. Stronger plumes will thus 
be found within the pulse, a fact corroborated by the larger 
kinetic energy densities found there (not shown). Both as- 
cending and descending plumes tilt leftward in such a right- 
traveling wave. In contrast to R T = 1.15/? ( /* (not shown) the 
solute contours possess clear regions of localized gradient 
inversion inside the pulses, indicating a local destabilization 
by the solutal field. The temperature also reveals regions of 
localized gradient inversion, but mild in comparison to those 
of the solute. No such regions of inverted gradients are 
present outside the pulse at this value of the Rayleigh num- 
ber. Outside the pulse kinetic energy is mostly dissipated in 
boundary layers near the top and bottom, where the no-slip 
boundary conditions force the velocity to vanish, while in the 
pulse region more kinetic energy is dissipated within the 
bulk of the fluid, hinting at a more active participation of 
internal friction mechanisms between the counterpropagating 
upward and downward plumes of fluid. 

Figure 4(a) shows the Nusselt number inside and outside 
the pulse. The values on the right correspond to the averaged 
Nusselt number N T and to the reduced Nusselt numbers in- 
side the pulse (N T p ) and in the background (N T b ). The 
Nusselt number N T (x) oscillates vigorously inside the pulse. 
When averaging these curves to obtain the reduced Nusselt 
number, a larger value is obtained inside the pulse region 
than outside, indicating that more efficient transport of heat 
takes place inside the pulse. Up to an overall scale the solutal 
Nusselt number N s (x) behaves in the same way and is omit- 
ted. Thus although the region occupied by the pulse is 19% 
of the total length of the convection domain, 24% of the total 
heat transport and 26% of the total solute transport is accom- 
plished within the pulse. 

2. R t =L 23R { t o) 

As the Rayleigh number is increased, more pulses are 
seen in mature states. Figure 2(b) shows a state with two 
right-traveling pulses at R r — 1.23 R { jK The number of roll 
pairs in the pattern remains at 26, the same as in the one- 
pulse case. The overall portion of the physical domain occu- 
pied by the pulses is now 36%, and the local wavelength of 
the rolls inside the pulses (\ r 1.21) is roughly 33% of that 

of the rolls outside the pulses (\ r ^3.68). The rolls move 
to the right at an average velocity of v r ^ 13.11, and both 
pulses travel to the right with the same velocity v p ^5A2, 
approximately 41% of the average velocity of the rolls. Al- 
though both pulses are of roughly the same length and re- 
semble the single pulse of the previous case, they are not 
equally spaced. The distance between the two pulses seems 
to be an equilibrium separation that is unaffected by small 
perturbations. 

The corresponding fields and Nusselt numbers are shown 
in Figs. 3(b) and 4(b). The two pulses now occupy 36% of 
the domain and carry 47% of the heat transport and 50% of 
the solute transport. The pulses are characterized by large 
regions of localized gradient inversion in the solutal field and 
milder ones in the temperature field. No such gradient inver- 
sion is observed outside them. 
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t = 30.0 f 


(top panels) and lateral heat flux density (bottom panels, continuous 
lines) and lateral solute flux density (bottom panels, dotted lines) 
for R T = 1.23 Rj >] . Left panels correspond to horizontal averages 
performed over the whole domain. Center (right) panels depict 
these same averages performed only over the pulse (background) 
region, as indicated by the reference panel at the bottom. 

piecewise constant, with relatively sharp transitions near the 
pulse boundaries. The rapid albeit smooth transition between 
them produces a clear resolution of the pulse boundaries 
when the time slices are taken together as in Fig. 6. This is so 
for the roll wavelengths as well [see Fig. 10(b)]. The wave- 


0 x —* 64 

FIG. 4. Thermal Nusselt numbers at a particular instant in time. 
Vertical lines indicate roll boundaries (dotted), and pulse bound- 
aries (solid). The bottom panel shows the pointwise behavior of 
N r (x) oscillating around the total thermal Nusselt number N r indi- 
cated at the right. In the top panel the numbers on the right indicate 
the total Nusselt number (middle), and the reduced Nusselt numbers 
inside (top) and outside (bottom) the pulse region. Solutal Nusselt 
numbers behave in a qualitatively similar way (see Table I for a 
summary), (a) R T = 1. 19 R { j ] . (b) R r = 1 .23 Rj >} . (c) R r 

= 1.27 R ( r o) . 

Figure 5 shows the instantaneous profiles of the mean 
lateral velocity u(zJ) and lateral heat flux density Q T izJ) at 
Rj~ 1 .23/?^ , with overbars denoting horizontal averages 
over either the whole domain or over the pulse and back- 
ground regions. The profiles of the lateral solute flux density 
Qs(z,t) are shown as dashed lines in the second set of pan- 
els. The left column contains averages computed along the 
whole physical domain, the middle column shows averages 
computed within the pulse, and the right column those com- 
puted in the background region ahead of the pulse. The con- 
tour plot of the stream function below the two panels iden- 
tifies the regions used in the center and right column 
averages. The mean lateral velocity u changes direction three 
times in the vertical, but the overall mass flux m is directed 
toward the left. 

Figure 6 shows the local phase velocity as a function of 
position in the domain in the form of a space-time plot. The 
pulses are associated with localized deficits in the phase ve- 
locity. Like the wavelength the phase velocity is nearly 
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FIG. 6. The roll phase velocity v r as a function of position in a 
two-pulse state at R r = 1.23 Rj-^ with 26 roll pairs at three times 
f = 25, 27, 29. Small solid dots indicate computed values. Vertical 
lines indicate roll boundaries (dotted), and pulse boundaries (solid). 
The reference stream function panels correspond to the first (bot- 
tom) and last (top) phase velocity plots. Values at the right side of 
each panel are (from top to bottom) the average roll velocity in the 
background, the roll velocity at the center of the pulse, and the 
pulse velocity. 
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lengths inside the pulses and outside are nearly constant, and 
both pulses show a clear difference between their leading 
and trailing ends. The presence of the second pulse is neces- 
sarily responsible for the larger wavelength of the rolls out- 
side the pulse region. Associated with the trailing interfaces 
one can now discern small superposed oscillations in both 
the local wavelength and phase velocity [cf. Fig. 10(b)]. The 
asymmetry between the leading and trailing ends of the 
pulses is present in all diagnostic quantities and is character- 
istic of all the pulses discussed here; it provides the key to 
the theory discussed in Sec. IV. 

3. R t =L27R^ 

Figure 2(c) also shows two rightward moving pulses, this 
time for R r = \ 21Rj }) . In this case 32 pairs of convecting 
rolls are present with an average wavelength of A. r **2.45. 
Although this state was obtained from the previous one 
(R r = \.23Rj y) ) by slowly increasing the thermal driving 
force, the total number of rolls in the pattern has changed. 
Therefore this state of the fluid motion is not a small pertur- 
bation of the previous one. The local wavelength of the rolls 
inside the pulses (\ rp ~\M) is still roughly 32% of that 
outside them (\ r ^3.43), but the pulse spacing is now uni- 
form. The average velocity of the rolls is v r ** 10.66, and that 
of the pulses is 17^4.49, approximately 42% of the average 
velocity of the rolls. 

The shape of the rolls in the central part of the pulses 
remains quite homogeneous, but the contour lines outside the 
pulses appear more disordered, suggestive of lateral erosion 
at the boundaries of the pulses where the stream function 
now attains larger values than inside. The trailing pulse is 
slightly larger than the leading one but the roll wavelengths 
away from the interfaces remain constant and similar in both 
pulses. As in the previous case, the relative pulse separation 
remains constant once a stationary regime is achieved. 

Figures 3(c) and 4(c) show that the trends observed at 
R t = \A9 Rj }) and R T = \.23Rj y) are enhanced. The pulses 
are still characterized by strong convection and solute gradi- 
ent inversion, but now regions of local gradient inversion in 
the solutal field can also be found outside the pulses. The 
tight packing of the convection rolls inside the pulses makes 
these gradient inversion regions less extended inside the 
pulses than outside, where the bent plumes can be ejected 
horizontally for larger distances, as can be seen in the solutal 
field 5 in Fig. 3(c). Inside the pulses, on the other hand, they 
are forced to recirculate and the gradient inversion is weaker. 
The stream function contours inside the pulses retain their 
symmetry; outside the pulses the maximum of the stream 
function is displaced laterally from the geometric center of 
the rolls, creating large horizontal gradients in 'T and there- 
fore large vertical velocities characteristic of strong ascend- 
ing and descending plumes of fluid. This accounts for the 
regions of gradient inversion in the solute outside the pulses, 
and is reflected in larger kinetic energy density there. The 
two pulses occupy 44% of the domain and transport 54% of 
the heat and 57% of the solute. 

4. R t =1.55R \ r o) 

For even higher values of the Rayleigh number pulses of 
constant form are no longer present. Figure 7 shows a trav- 



0 X ► 64 


FIG. 7. As for Fig. 2, but showing a traveling wave at R r 
- 1.55 Rj? } ■ Two nearly stationary pulses are present. Spatiotem- 
poral defects form both inside the pulses and in the background; 
phase is no longer conserved. 

eling wave pattern at R T = 1.5 5Rj >] . The pulses are replaced 
by a couple of almost stationary regions of very small hori- 
zontal extent that the pulses have to traverse. These “barri- 
ers” serve as phase sinks: two events leading to the disap- 
pearance of a roll pair can be discerned, one in each 
“barrier.” In contrast the background region serves as a 
phase source: two events creating an additional roll pair can 
be seen in the time interval shown. It appears that the num- 
ber of roll pairs remains statistically stationary, with the 
number of roll-splitting events in the background balanced 
by roll-annihilation events in the slow-moving “barriers.” 
The presence of these spatiotemporal defects in the pattern 
indicates that phase is no longer conserved. 

5. Mutated states 

We have found that it is possible to generate multiple, 
numerically stable states for fixed parameter values by a pro- 
cess we call mutation. In this process selected parts of a 
stable wave train were compressed while others were ex- 
panded, with the overall number of roll pairs kept un- 
changed. Such an initial condition obtained by “gluing” the 
pieces back together was then evolved by time stepping for- 
ward in time. If it was the pulse that was compressed the 
number of rolls in the pulse decreased rapidly as the system 
attempted to recover its equilibrium pulse wavelength. In 
contrast, when the background rolls were compressed a new 
pulse would often develop, again on a fast time scale. In all 
cases the resulting state was evolved until a numerically 
stable final state was reached. In this way different stable 
states could be manufactured for the same parameter values. 
The new states would typically differ from the original state 
in the overall roll number, the separation of the pulses, and in 
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FIG. 8. As for Fig. 2, but showing for R r ~ 1.23/?^ 1 (a) a three- 
pulse state, and (b) a two-pulse state, both obtained by mutation of 
the two-pulse state shown in Fig. 2(b) for the same R T . The two- 
pulse state differs from the 26 roll pair state shown there, although 
it has the same number of rolls. 


the number of oscillations in their tails. They may even differ 
in the overall number of pulses. We show two such states at 
R t = \.23Rj }) in Fig. 8 and summarize their properties in 
Table I. Other stable states have been found by starting from 
a stable state and either increasing or decreasing R r . For 
example, at R r = 1.18/?-/'* we found a single pulse using the 
numerically stable pattern found for R r = \A9Rj >} as initial 


configuration, but three pulses when starting from the nu- 
merically stable pattern found for /? r = 1.15/? ( / ) (not shown). 

6. Summary of transport properties 

The transport properties of the pulses are summarized in 
Tables I— III. In all cases the total thermal and solutal Nusselt 
numbers N r and N s are increasing functions of R T . Within 
the states sharing the same Rayleigh number ( R T 
= 1.23 /? ( / ;) ), the Nusselt numbers increase with the number 
of rolls in the state. The same holds for the reduced Nusselt 
numbers inside the pulses, N Tp and N S p , but not for the 
reduced Nusselt numbers in the background (see Table I), 
indicating a gradual increase in convection efficiency inside 
the pulses but not in the background region. 

For R t = 1.23/?j? ) the two states with 26 roll pairs and two 
pulses share almost identical transport properties, being dif- 
ferentiated only by the relative size of the pulses and their 
mutual separation. Despite these differences the total extent 
of the pulse and the background regions is roughly the same 
in both cases, with the pulse widths comparable in the first, 
and the two background regions comparable in the second. 

When phase is conserved (as for R T = 1 . 1 9R { f ) , 1.23 R { f\ 
\.21Rj >y ), the pulses are very active dynamical structures 
that play a substantial role in convective transport. In this 
regime the proportion of the kinetic energy, of the enstrophy, 
and also of the vertical heat and solute transports attributed 
to the pulses exceeds the relative portion of the domain oc- 
cupied by the pulses. Thus pulses enhance the overall trans- 
port properties of the system (cf. Table II). Moreover, when 
two pulses are present, the fraction of the domain they oc- 
cupy increases both with R r and with the number of rolls in 
the pattern. In contrast, when phase is not conserved (as for 
R 7 = 1.31 R\ 0) and 1.35/? ( r 0> ) the proportional heat and solute 
transports by the pulses are comparable to the proportional 
extent they occupy in the horizontal, indicating that the pres- 
ence of the pulses no longer increases transport efficiency. 

In all cases studied the average velocity of the rolls (i/ r ) 
decreases with increasing Rayleigh number while the overall 
average roll length (X,.) increases. The roll velocity within a 
pulse (v r p ) also decreases with increasing Rayleigh number, 
while in states corresponding to the same Rayleigh number 
(Rj= 1.23 Rj^) it also decreases with increasing number of 
rolls (cf. Table IV). The velocity of the pulses is much less 
than the phase velocity of the long-wavelength modulation 
seen at small Rayleigh numbers (R T = \A5R ( j ] ) y indicating 
that the pulses are fully nonlinear states. 

Table III compares the lateral fluxes of mass, heat, and 
solute obtained by averaging the horizontal fluxes corre- 
sponding to 50 snapshots spanning five units of time. There 
is a small lateral mass flux m, which for the Rayleigh num- 
ber values analyzed ranges between -0.01 and -0.06, 
namely, opposite to the motion of the rolls. This lateral mass 
flux originates from an asymmetry and tilting of the rolls. 
The lateral transports of heat and solute are also directed 
toward the left in these cases. These lateral fluxes are small 
when compared to the vertical ones (as judged by the ratios 
Qj!N t or Q s fN s ), namely, between 60 and 120 times 
smaller for the heat flux, and between 25 and 60 times 
smaller for the solute flux (cf. Tables I and III). In all cases 
the lateral transports show large relative deviations and their 
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TABLE I. The thermal and solutal Nusselt numbers computed for the whole domain (N T ,N S ), within the pulses (N Tp ,N s p ), and in the 
background (N Ttb ,A f s ^ b ). Numbers in parentheses indicate standard deviations. 



Roll pairs 

Pulses 

N t 

N T . n 

A T,b 

Ns 

Ns, P 

Ns.b 

1.15 

26 

0 

1.43 



1.52 






(9.5XKT 5 ) 



( 1.4X 10“ 4 ) 



1.19 

26 

1 

1.56 

1.99 

1.46 

1.71 

2.32 

1.56 
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(3.5X I0“ 2 ) 

(1.0X 10“ 3 ) 

(5.5x 10“ ') 

(1.1X10“') 

1.23 

22 a 

2 

1.57 

2.06 

1.39 

1.73 

2.42 

1.48 




( 2.1 x io~ 3 ) 

(1.9X 10“') 

(6.5X10 -2 ) 

(4.2X10 3 ) 

(5.9X10 ‘) 

(2.0X 10“') 


24“ 

2 

1.62 

2.13 

1.39 

1.81 

2.54 

1.48 




(2.7X 10' 3 ) 

(6.5X10 -2 ) 

(2.8X 10“ 2 ) 

(5.6X 10 3 ) 

(2.1 X 10 1 ) 

(8.7X 10 2 ) 


26 

2 

1.67 

2.18 

1.38 

1.89 

2.63 

1.47 




(3.4X 10“ 3 ) 

(1.2X 10“') 

(6.6X10 -2 ) 

(7.5X 10“ 3 ) 

(3.7X10“') 

(2.1 X 10“') 


26 ab 

2 

1.67 

2.17 

1.38 

1.89 

2.62 

1.47 




(2.6X 10“ 3 ) 

(8.9X 10“ 2 ) 

(4.6X 10 2 ) 

(6.1 x 10 3 ) 

(2.8X10“') 

( 1.4X 10“') 


33 a 

3 

1.88 

2.24 

1.44 

2.20 

2.75 

1.53 
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( 1.2X 10“ 1 ) 
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(3.7X10“') 
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1.27 

32 

2 

1.95 

2.40 
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2.34 

3.05 
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2.73 
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Corresponds to mutated states. 







^This two-pulse, 26 roll pair state differs from the one shown in Fig. 2(b) in the relative size of each of the pulses and 

in the interpuls 


separation. 

c Phase is not conserved for these values of R T . The number of roll pairs fluctuates or is not precisely defined; the value reported is the most 

frequent in the time interval analyzed. 

d Pulses are no longer seen for this value ot R T . 


dependence on the Rayleigh number remains unclear. 

B. Anatomy of a pulse 

Thus far we have characterized our solutions either by 
means of snapshots at particular instants ot time or by means 
of horizontal averages. We now describe the properties of 
comoving time averages of temperature, solute, stream func- 
tion, and vorticity, as well as of the thermal and total buoy- 
ancy forces defined by 

A T -«r0. f B S) =~ R sZ. /fl s // , +// ) =Kre-*s2- 

(17) 

Time averaging of spatiotemporal patterns has proved to be a 
successful tool in analyzing complicated pattern- forming 
systems. In some cases it can reveal regular structures under- 
lying disordered ones, even well into spatiotemporal ly cha- 
otic regimes [17]. Indeed Barten et al. [16,7] have used time 
averages of fields to deduce the role played by the concen- 
tration field in holding together two different coexisting do- 
mains. 

Time averages can be performed in the laboratory frame 
or from the reference point of a moving observer. The former 
provide essentially the same information as the instantaneous 
horizontal averages of the preceding section, and will not be 
considered further. We denote the latter by ( )„, with the 


subscript u indicating the velocity at which the observer is 
moving. Carefully chosen time averages filter out small-scale 
features of the fields and reveal their long-range components. 
In particular, the time averages recorded by an observer co- 
moving with the pulses allows one to infer what the pulses 
are actually “seeing’ 1 when they drift along the domain. 

A series of 50 to 100 consecutive snapshots of the fields 
was taken at equal time intervals of 0.1, spanning from five 
to ten units of time. These fields were suitably shifted and 
then averaged together. Thus, if F(jc,z, 0 represents any of 
the fields, then 

1 i « 

<F(jc, 2)> b =- I F(x-vt,z,t)dt~- X F(x-vtj,z,tj). 

T Jo n j = 1 

08 ) 

The resulting structures were further smoothed with a ta- 
pered low-pass filter. 

Figure 9 shows the time-averaged fields at /? 7 = 1.23F ( / ) 
(with 26 roll pairs in the pattern) as seen by a comoving 
observer traveling at the pulse velocity v p ~5A2. The aver- 
age reveals clearly the pulse location. The temperature and 
solute time averages possess smooth gradient variations, 
from T=S= 1 at the bottom, to T=S = Q on top. Inside the 
pulse region, the vertical gradients of temperature and solute 
are smaller at midlayer than in the surrounding fluid. Those 
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TABLE II. The ratio A p /T denotes the fraction of the domain occupied by the pulses. Shown also are the 
fractions of the kinetic energy E p IE and enstrophy Z p !Z within the pulse regions relative to those quantities 
or the whole domain, and the proportion of the vertical heat flux carried by the pulses compared to the total, 

namely, N T p A p /N T T, and the analogous proportion of the vertical solute flux. Numbers in parentheses 
indicate standard deviations. 


r t / r ( t ° } 

Roll pairs 

Pulses 

Vr 

E p /E 

z p /z 

N Tp A p m T r 

N s _ p A„IN s T 

1.19 

26 
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0.19 
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0.54 
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J Corresponds to mutated states. 






This two-pulse, 26 roll pair state differs from the one shown 

in Fig. 2(b) in the relative size of each of the 


pulses and in the interpulse separation. 

Phase is not conserved for these values of R T . The number of roll pairs fluctuates or is not precisely defined; 
the value reported is the most frequent in the time interval analyzed. 


same gradients are, however, very large close to the top and 
bottom boundaries, accounting for the larger reduced Nusselt 
numbers inside the pulses. For ease of comparison the total 
and thermal buoyancy force fields, as in Eq. (17), are plotted 
using the same set of isobuoyancy contours. Thermal and 
solutal (not shown) buoyancy forces possess rather feature- 
less profiles, much like the temperature and solute averages, 
except for a couple of pronounced maxima and minima that 
reveal the location of the pulse. Both combine to give a 
complicated structure to the total buoyancy force in the pulse 
region. 

The stream function panel in Fig. 9 shows an average 
long-range flow that circulates along the whole domain from 
right to left, with streamlines wiggling up and down through 
the bulk of the fluid. This is a familiar property of traveling 
wave streamlines in the comoving frame [18]. In the labora- 
tory frame the magnitude of the lateral mean flow spanning 
the domain was too weak to reveal itself readily in the stream 
function contour plots [cf. Fig. 3(b)]. The striking feature of 
this averaged stream function is the presence of groups of 
vortices signaling the position of the pulses. There are two 
pairs of vortices for each pulse, one at the leading edge of the 
pulse and the other at the trailing edge, each pair with one 
vortex above and one below midlevel in the domain. The 
vortex pair at the trailing edge (on the left) is the strongest of 
the two, with the one at the top rotating counterclockwise 
and the lower one clockwise, both embraced by the flow 
lines of the long-range flow. The pair at the leading edge is 
weaker and covers a larger extent of the pulse area, with a 


clockwise vortex on top and a counterclockwise one at the 
bottom. 

It is not surprising to find a cluster of vortices traveling 
with the pulses. We should recall that in all cases more effi- 
cient transport was achieved in the pulse region; the back- 
ground might have stronger plumes, but the ones inside the 
pulse are closer together. The plumes provide a coherent 
means for transporting both heat and solute. These transports 
depend on the correlation between the vertical component of 
the velocity vc and the temperature and solutal fields T and 
S. The pulses are therefore regions of enhanced correlation 
between these fields and also of enhanced coherence among 
the plumes. The resulting spatial nonuniformity implies the 
presence of large-scale horizontal gradients of w and hence 
the presence of large-scale vorticity. These vortices cannot 
be seen in Fig. 2 because they are much weaker than the 
convection rolls. 

IV. PHASE EQUATION APPROACH TO THE PULSE 
FORMATION PROBLEM 

The pulses described in Sec. Ill were found to be very 
stable against random noise perturbations added to the fields 
during the numerical simulation. When two pulses were 
present, their separation and propagation speeds were unaf- 
fected by this type of perturbation, as long as the noise in- 
tensity remained below certain threshold values. When the 
noise level is set above the threshold the pattern is destroyed 
and the process of roll formation resumes almost from 
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TABLE III. Time-averaged lateral mass flux m, the lateral heat flux Q T which represents the vertical 
integral of Q^zh and similarly the lateral solute flux Q s . Numbers in parentheses indicate standard devia- 


tions. 
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a Corresponds to mutated states. 

^This two-pulse, 26 roll pair state differs from the one shown in Fig. 2(b) in the relative size ot each o the 


pulses and in the interpulse separation. 

c Phase is not conserved for these values of R T . The number of roll pairs fluctuates or is not precisely defined; 
the value reported is the most frequent in the time interval analyzed. 
d Pulses are no longer seen for this value of Rj. 


scratch, without memory of its past state. In contrast, when 
an existing state is mutated in the manner described in Sec. 
Ill A 5 the system usually relaxes into a new stable state. The 
numerical evidence thus favors the coexistence of multiple 
stable states. In this section we seek to elucidate (a) the na- 
ture of an individual pulse (its shape, propagation speed, 
etc.), (b) the origin of the double-pulse states described in 
Sec. Ill, and (c) the apparent multiplicity of pulse states al- 
luded to above. 


A. Kinematics 

In all three cases discussed in detail in Sec. Ill the phase 
of the pattern was found to be conserved, that is, the number 
of convecting rolls stayed constant through time, without 
merging or splitting of rolls. Our theoretical approach takes 
this observation as its starting point. We think of the travel- 
ing wave in the form ^(kx— <p(x,t),z)-, where <p(x,t) 

is the phase of the wave. For a uniform wave train the phase 
is constant. However, a spatially dependent phase provides a 
convenient description of nonuniform wave trains with wave 
number q(xj) and frequency fl(x,r) that depend on posi- 
tion. 


The simulations indicate a uniformly moving pulse, moving 
with speed v p . The modulation in wave number and fre- 
quency thus takes the form 

(p(xj) = (p(x-v p t). (20) 

Thus 

q = k + <p\ fl= w+ Vpip', (21) 

where the prime denotes differentiation with respect to argu- 
ment, i.e., <p' = d(pld£, £=x~v p t. Since the simulations in- 
dicate that (p(xj) is defined for all (x y t) we have the kine- 
matic relation 


dq d£l 
dt + dx 


= 0 


( 22 ) 


that expresses conservation of waves [19]. In the frame of the 
pulse this relation takes the form 



dq dfl 


(23) 


dip 


q=k+ ^’ 


dip 

~dt 


U\q-k\<\, | ft -a> |^1 while the pulse interfaces (leading 
(19 ^ or trailing) have 0( 1 ) thickness the appearance of the pulse 
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TABLE IV. Comparison of computed pulse velocity v p with values predicted from the jump condition 
(25). See Sec. Ill A for symbol definition. When phase is conserved, typical measurement uncertainties for 
wavelengths and phase velocities are ±0.01, and for pulse velocities ±0.02. 


Rj/R^ 

Roll pairs 

Pulses 

K 

X r .6 


Vr 

V r,b 

V r.P 


v (p)b 

v p 

1.15 

26 

0 

2.00 



15.15 



1 1.65 c 


1.19 

26 

1 

2.12 

2.86 

1.31 

14.21 

15.01 

10.05 

5.88 

5.86 

1.23 

22 d 

2 

2.20 

3.72 

1.42 

14.03 

15.27 

9.33 

5.66 

5.66 


24 d 

2 

2.28 

3.65 

1.26 

13.46 

15.02 

8.69 

5.37 

5.38 


26 

2 

2.38 

3.68 

1.21 

13.11 

15.05 

8.57 

5.42 

5.41 


26 d,e 

2 

2.38 

3.71 

1.19 

13.11 

15.07 

8.52 

5.44 

5.43 


33 d 

3 

2.44 

3.65 

1.20 

11.66 

15.11 

8.46 

5.39 

5.21 

1.27 

32 

2 

2.45 

3.43 

LI l 

10.66 

12.83 

7.24 

4.49 

4.57 

1.31 

46 f 

2 

2.27 

3.33 

0.99 

7.55 

11.77 

5.99 

3.63 

3.54 

1.35 

45 f 

2 

2.32 

3.03 

0.92 

6.97 

10.06 

5.14 

3.29 

3.01 


d Measured pulse velocity values. 

Predicted pulse velocity values. 

Corresponds to small-amplitude long-wavelength modulations. 

Corresponds to mutated states. 

e This two-pulse 26 roll pair state differs from the one shown in Fig. 2(b) in the relative size of each of the 
pulses and in the interpulse separation. 

f Phase is not conserved for these values of R T . The number of roll pairs fluctuates or is not precisely defined; 
the value reported is the most frequent in the time interval analyzed. 


is very abrupt and can then be viewed as a shock. Suppose 
that such a shock (again leading or trailing) forms at £= £ 0 in 
the pulse frame. Integrating across the shock and noting that 
the shock is, by construction, stationary yields the “jump” 
condition 


R t = 1.23 R[ o) 



0 X ► 64 


FIG. 9. Time-averaged fields in a frame traveling to the right 
with the pulse velocity for R T = 1.23 Rj ] . Shown are the averaged 
temperature (T) Vf> , stream function (¥)„ , total buoyancy force 
( /b)v p * anc * thermal contribution to the total buoyancy force 
(fjp) Vp * obtained from 100 snapshots spanning ten units of time. 
The beginning 'P field is included to specify the change in the 
vertical scale. The corresponding time average of the solutal con- 
centration S is similar to that of the temperature T and has been 
omitted. The solutal contribution to the total buoyancy force 
</b S \ closely resembles -{fP) v and has also been omitted. 


v p &q = A£l, (24) 

where A q — q + ~q~, A fl S5 fl + — fl“ represent the jumps in 
the wave number and frequency across the shock. We can 
rewrite this relation in the more convenient form 

k + v~ — k~v + 

"' = “V-k- • < 25) 

where v = Cl/q is the local phase velocity and k = 2ir/q the 
local wavelength. The predictions from this formula are 
compared with measured values of v p in Table IV. Since the 
pulses do not move completely rigidly [their widths mea- 
sured by tracking regions where \,(;c)<\ r exhibit small os- 
cillations due to the passage of rolls in and out of the pulse] 
the values of v p computed from Eq. (24) are not completely 
time independent. Consequently the values of k r b and v r b 
listed in Table IV represent the average of their instanta- 
neous values, measured frame by frame. Similarly, the mea- 
sured pulse velocity v p also represents a time average. None- 
theless the average of the instantaneous values of v p 
computed from the jump condition (25) using the instanta- 
neous values of k pb and v pb shows remarkable agreement 
with the pulse velocity measured directly. 

B. Dynamics 

The phase <p of the wave V(kx - + (p(xj),z) satisfies 
an evolution equation called a phase equation. This equation 
follows from a gradient expansion (d x <p<\, 1 ) and the 

symmetries of the underlying wave train. In the present case 
Galilean symmetry is absent (because of the no-slip bound- 
ary conditions imposed at z = 0,l); consequently the wave 
train has only two symmetries, a continuous translation sym- 
metry, 


<p— + const. 


(26) 
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indicating the fact that the absolute value of the phase is 
arbitrary, and the discrete “mirror-glide” symmetry [2,14— 
16]. If the wave train has period litlk this symmetry acts by 

'P— + — 4L <£—►<£+— . (27) 


Consequently the mirror-glide symmetry has no effect on the 
structure of the phase equation. The most general evolution 
equation for the phase is therefore 

dip dip d 2 if d*ip d 4 ip 

dt +C * dx~ U ~dx I ^ "dir C \ dx / 


— 


dip d~<p 
dx dx 2 




/ d 2 ip \ 2 dip d* ip I d<p\ 4 

h Mr ~k + *** * 

\ dx~ / dx dx ' \ dx I 

(28) 


correct to fourth order in spatial derivatives, cf. [20]. Here c g 
is the group velocity of the perturbations at onset in the 
laboratory frame. This equation typically exhibits solutions 
in the form of spatiotemporal chaos, called phase turbulence 
[21]. Such spatiotemporal chaos occurs when 0<a<^l, and 
is to be distinguished from the appearance ol “phase slips 
at which the phase description breaks down. These typically 
occur for larger values of a. As already mentioned, the pres- 
ence of phase slips invalidates the phase description; this is 
not so when phase turbulence is present. In order to model 
the parameter regime R r — ( \ A9- \ .21)Rj >] in which stable 
pulses are observed, and phase turbulence absent, we restrict 
the range of possible behavior by imposing the requirement 
that the right side of Eq. (28) be in conservation form, i.e., 
</=/=/ = 0, h = k. Thus 


dip dip dj 

dt + C H dx dx ’ 


(29) 


where 

j = aip x +b<p xx + c<p XXK + e<p 1 x + g<p\ + hip x ip xx +-" . 

(30) 


tional and consequently admits only time-independent states 
in the long time limit. This is also the case in the closely 
related Chapman-Proctor equation ( c g = b = e = h = 0 ) [25]. 
With appropriate boundary conditions all these equations 
have steady pulselike solutions, as discussed in detail by 
Deissler, Lee, and Brand [23]. However, because they have 
variational structure all stable pulses are necessarily subcriti - 
cal. This is because the equilibrium connected by the pulse 
must itself be stable. As a result the pulse formation process 
is hysteretic. In nonvariational systems this is not necessarily 
the case. Indeed pulses have been observed at supercritical 
parameter values in experiments on binary fluid convection 
[26]. Moreover Schopf and Kramer [27] found stable pulses 
in the complex Ginzburg-Landau equation with a destabiliz- 
ing cubic nonlinearity at supercritical values of the bifurca- 
tion parameter, while in systems with a nonzero group ve- 
locity pulses can also be stable in the supercritical regime, 
provided the basic state is only convectively, but not abso- 
lutely, unstable. This is the regime of interest in the present 
connection. 

From Eq. (29) we see immediately that infinitesimal 
phase perturbations obey the dispersion relation 

s = aQ 2 -iQ 4 + --+iQ(c >t -bQ 2 +--), (31) 

where ip(xj)~exp(st-iQx). Thus phase perturbations grow 
when aQ 2 ~cQ 4 + • • • >0. Since the anticipated instability is 
of long wavelength (Q < ^\ ) we suppose that a and c are both 
positive. The observed pulses are the outg rowt h of the result- 
ing convective instability (0<a<2. \2\tbc g if bc g > 0 , 0 
<a< cc if bc g <0 , when c = 0). We have seen that the insta- 
bility evolves into a fully nonlinear pulse moving with speed 
v p . Even though the phase equation (29) was obtained by a 
gradient expansion we shall find it to be useful as a model 
even for the fully nonlinear pulses for which Q = 0( 1 ) seen 
in the simulations. A steadily moving pulse moving with 
speed v p to the right satisfies the equation 

d{ [{c H -v p )ip + aip£+bip^+Cip^+eip 2 £ + g<p] + hip{ip^ 

-\ ] = 0. (32) 


Thus 


The motivation for studying Eq. (29) instead of Eq. (28) is 
threefold. First of all, it allows us to make contact with re- 
lated work on confined states arising through phase instabili- 
ties of periodic patterns of steady rolls [22,23]. Second, it 
allows us to reduce the number of coefficients in the descrip- 
tion. Finally, the phase equation in the pulse frame can be 
integrated once, yielding a third-order equation for the pulse 
shape. Given the fact that the observed pulses are character- 
ized by a monotonic rise and an oscillatory trailing end such 
a description is the minimal one consistent with the observa- 
tions. Moreover, as discussed further below', the observed 
pulse shapes are completely consistent with such a third- 
order description. 

Several special cases of Eq. (29) have already been stud- 
ied. For steady rolls the corresponding equation (c g = b = h 
= ()) is of B- E type in the phase equation classification of 
Kuramoto [24]. The resulting equation has a Liapunov func- 


vtp + atpt+b(p£i+c<PK£+e<p\ + gip\+h<pt<ptf+ 0 , 

' (33) 

where v = c g -v p and the constant of integration has been set 
to zero. In the following we think of Eq. (33) as a nonlinear 
eigenvalue problem for the speed v, and hence for the pulse 
speed v p . As discussed further below the formulation of the 
boundary conditions for this problem depends on whether it 
is posed on the whole real line or on a finite (periodic) do- 
main. If, in addition, the number of roll pairs in the perturbed 
and unperturbed wave trains is unchanged as a result of the 
formation of the pulse, there is no net phase jump A ip across 
the domain; this is not so if the number of roll pairs differs. 
The remaining boundary conditions for Eq. (33) are either 
that both <p£ and (p H vanish as ± 30 (unbounded domain) 
or that both are periodic with period T (bounded domain). 
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C. Pulses as homoclinic orbits 

We begin with the eigenvalue problem (33) on the real 
line, and use it to draw a number of general conclusions 
about the pulse shape and other properties. To this end we 
set c= 1 and rewrite Eq. (33) as a dynamical system in £. We 
write * for <p, so that Eq. (33) becomes 

.v' = Z> (34) 

z'=-vx- ay ~bz~ ey 2 ~ g y 3 - hyz, 

with the prime denoting dldg. This dynamical system de- 
scribes the spatial structure of the wave number modulation 
which is given by y({j) in the present notation. It is this 
quantity that can be compared with our simulations. It is a 
dissipative system, 


dx' 

~f- 

dx 


dy r>Z 


— b — hy. 


(35) 


does so as a decaying exponential without superposed oscil- 
lations. Thus the pulse shape is intrinsically asymmetric: 
when v p <c g the leading edge of the pulse is monotonic 
while its trailing edge is an exponentially decreasing oscilla- 
tion. The opposite is the case when v p >c g , in which case 
the pulse leads with an oscillatory precursor. In case (ii) 
which applies for larger |y|, e.g., for v p <c g (v p »>c g ), the 
oscillatory tail (precursor) is absent, with the interface re- 
placed by a monotonic one. The pulse remains generically 
asymmetric, however. These conclusions continue to hold 
when 0<|6|<^ I and 

4(i 2 b 


s 2.3 =/*' 


5/J.-+ v 2 b 
9 /z 2 + v 2 2 


, 3/u 2 — v 2 (jib 


(38) 


with a single equilibrium at the origin corresponding to the 
unperturbed wave train. 

Equation (34) is in the form of a particular codimension- 
three bifurcation problem [28], When v -a = b= 0 the fixed 
point (0,0,0) of Eq. (34) has three zero eigenvalues; the pa- 
rameters v,a,b can therefore be thought of as unfolding pa- 
rameters that “unfold” this triple zero degeneracy. The 
pulses that are of interest correspond to homoclinic orbits of 
Eq. (34) that start from the origin at £= — oo and return to it 
at £= -foe. The excursion of such an orbit (in £ space) away 
from the origin defines the shape of the pulse. Evidently, the 
question of whether such pulses exist is a global one. How- 
ever, a necessary condition for the existence of a pulse is that 
(0,0,0) is a generalized saddle in £, with at least one unstable 
eigenvector and at least one stable one. 

The linearization of Eq. (33) around the origin yields the 
characteristic equation 

s* + bs 2 + as + v = 0, (36) 

where, for phase instability, a>0. Thus if Res>0 the solu- 
tion grows with increasing £ but it decays with £ if Res<0. 
To construct a pulse both types of solution are necessary. 
Thus when v>0 we need either b< 0 or ab - v <0 (or both) 
and when KOwe need either 6>0 or ab~v> 0 (or both) 
The simplest situation arises when the phase modulations are 
only weakly dispersive so that 6 = 0. The sum of the eigen- 
values 5|,j 2 ,-r3 is then zero, and hence either (i) *,= 
-2(i, s 2i3 =(i±iv, or (ii) .y 2 =-/z 2 , .v 3 = /z , 

+ (i 2 ■ In the former case 

a=v 2 - 3 M 2 >0, 2 (i((i 2 +v 2 ) = v. (37) 

This case applies when v and hence (i are not too large, and 
is the case of greatest interest to us. Note, in particular, that 
sgn/i=sgnu. Thus, when u>0, the origin has a one- 
dimensional stable manifold and a two-dimensional unstable 
one; the opposite is the case when u<0. Hence when v>0 
(i.e., v p <c g ) the trajectory leaves (0,0,0) at £=-oo as an 
exponentially growing oscillation; when it returns as £— >oo it 


The numerical simulations provide valuable information 
about the pulse shape, and in particular about the wave num- 
ber and frequency modulation due to the presence of a pulse. 
In Fig. 10 we show the wavelengths and phase velocities of 
the convection rolls at R T = 1.19/?'/’, R T =\.23R { P \ and 
R t = 1 ,21R { j ) at one instant in time. Observe that the wave- 
length and phase velocity curves have similar appearance. 
Both inside the pulses and outside the phase velocity and 
wavelength distributions are piecewise constant, with fairly 
rapid transitions between two sets of values, and both curves 
resemble one another. For pulses on the real line the phase 
velocity is defined as v=Cl(x,t)/q(x,t) while the wave- 
length \ = 2 tt! q(x ,t) . For small <p' it then follows that 


V-V p = (v (39) 

*0 

where v 0 and X 0 are the phase velocity and wavelength of 
the background uniform wave train, i.e., v 0 =a>/k, \ 0 
= 277/6. Thus the phase velocity and wavelength are related 
linearly when phase gradients are small. Since u 0 —v p > 0 
(see Table IV) Eq. (39) indicates that the two curves should 
match without reflection. Figure 10 supports these conclu- 
sions; departures from the linear relation visible with in- 
creasing R t are not surprising since the phase gradients in- 
crease with increasing /? r . 

The comparison of Fig. 10 with the long- wavelength 
theory thus clearly favors our case (i) with v>0. The single 
pulse present at R T = L\9R { T o) leads with a monotonic front 
while its trailing edge decays exponentially at a different 
rate, with small superposed oscillations clearly visible. In the 
case of the two pulses found at R r = 1.23 R^ } both pulses 
have the same structure, with the oscillatory tail of the lead- 
ing pulse decaying towards negative f before turning into the 
leading front of the trailing pulse. At R T = 1.27 these 
oscillations are even more prominent, suggesting that jx is 
decreasing with R T and hence that the pulse velocity be- 
comes closer to the linear group velocity. 
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(a) R t = 119 R\ o) 



(b) R t = 1.23 R\ o) 



(c) R t = 1.27 Rt o) 



FIG. 10. Profiles of roll velocities and wavelengths as a function 
of position in the domain at one instant in time. The profiles travel 
almost rigidly to the right with speed u r ■ Vertical lines indicate roll 
boundaries (dotted), and pulse boundaries (solid). The dots corre- 
spond to computed values of wavelengths and phase velocities. Val- 
ues at the right of the phase velocity panels are, from top to bottom, 
the phase velocities in the background and inside the pulse, and 
pulse velocity v p . Values at the right of the wavelength panels are 
(top) the roll wavelengths outside the pulse and (bottom) inside the 
pulse, (a) R t = 1.19 . (b) R T = 1.23 R { ? . (c) R T = 1-27 R'j’ ■ 

However, despite this superficial agreement with the long- 
wave theory, the relation between the theoretical description 
of a pulse on the real line and a train of periodic pulses such 
as obtained in our numerical simulations is not straightfor- 
ward. This is because in a finite domain with periodic bound- 
ary conditions the pulse is a fully nonlinear state everywhere: 
even in the wings of the pulse the wave number modulation 
about the putative background state does not vanish. Thus 
the use of periodic boundary conditions has important impli- 
cations for the correct definition of the background state. To 
see this note that on the real line the definition of cp(.r,f) as 
the perturbation to the phase of a spatially uniform wave 
train and its localized nature imply that d£— <p(°°) 

- <p (-x) = Q. Thus the area under the wave number pertur- 


bation must be zero. On a finite but periodic domain the 
observed conservation of phase during the evolution of the 
phase instability likewise implies that J 0 <P^ d£= <f>(T) 
-<p(()) = 0. As can be seen from Fig. 10 the condition 
d£=0 defines a background (the underlying uniform 
wave train) that differs from the background wave number 
k h defined by the condition that |<p^- 1 > n l ^e largest 
possible fraction of the domain 0=s£=sr. Relative to this 
background Xq(<p^ — k b )d^—A>0, and hence k h — AIT 
<0. The wave number shift k h is finite because the phase 
perturbation is distributed over a finite domain and not the 
whole line; k b vanishes as T— Relative to the shifted 
background a pulse represents a net phase change, measured 
by A. In the following we therefore relate the observed 
pulses to solutions of Eq. (33) on with the bound- 

ary conditions <p(0) = 0, <p(0=A, and <p ff both periodic, 
with A determined from Fig. 10. Similar issues were encoun- 
tered already by Deissler, Lee, and brand [23] in their study 
of confined states described by phase equations with vana- 
tional structure. 


D. Determination of the dispersion relation 

In this section we use the shape of the pulse wings in Fig. 
10 to determine the coefficients ale, blc , and vie in the 
dispersion relation (31) on the assumption that the observed 
pulses approximate the homoclinic pulses described in the 
preceding section. In view of the finite area A and the 0(1) 
wavelength modulation in the observed pulses it is necessary 
that we first redefine the background wave number as dis- 
cussed above. To determine the dispersion relation govern- 
ing small wave number perturbations about this background 
we fit each interpulse region to a curve of the form 

f(x) = a 0 J ra l x + O2 e> ' ii ^~ e X2jr [a:iCos( or*) + fc^sin( trjx)]. 


connecting the leading edge ot one pulse to the trailing end 
of the preceding one. In cases where multiple pulses were 
present we focused on the largest interpulse region. We 
chose this approach instead of studying each pulse wing in- 
dividually because the characteristics of adjacent pulses typi- 
cally overlap. Even in the one-pulse case, where one would 
expect that the extent of the background region would make 
the leading and trailing pulse ends disconnected, the wig- 
gling of the trailing end propagates all the way up to the 
leading end. In all cases the fitted value of a 0 was very close 
to the measured wavelength of the rolls in the background. 
The parameter a { turned out to be zero in all the fits except 
in the one-pulse case ( R r = 1.19/f'j ’ ) in which a slight slope 
of the background is present, representing a wave number 
gradient between the pulse and its periodic image. The fitted 
values of X , , X 2 , and <r 2 provide estimates of the quantities 
s, , Re.s 2 , and Itm 2 . respectively, and hence of the roots a, , 
a 2 , F 2 of the characteristic equation (36). Table V lists these 
estimates. Such estimates improve with the interpulse sepa- 
ration. In a successful fit X | <0, X 2 >0. From these eigenval- 
ues we can obtain values of a, b, and v. To check on the 
consistency of these results we use the observed pulse speed 
v n to deduce from v the group velocity c K . Table VI lists the 
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TABLE V. Exponential fits to the largest interpulse intervals, 
computed by fitting the function (40) to the measured data by least 
squares techniques. Numbers in parentheses indicate standard de- 
viations. 


r t /r </> 

S] 

Res 2 

Irm 2 

1.19 

-0.82 

1.36 

0.90 


(0.08) 

(0.43) 

(0.61) 

1.23 

-0.82 

1.52 

1.05 


(0.45) 

(0.38) 

(0.46) 

1.23 a 

-0.72 

1.75 

1.12 


(0.11) 

(0.51) 

(0.41) 

1.23 b 

-0.92 

1.63 

0.74 


(0.14) 

(0.29) 

(0.48) 

I.23 c 

- 1.02 

1.51 

1.18 


(0.19) 

(0.65) 

(1.28) 

1.27 

-0.75 

1.17 

0.76 


(0.09) 

(0.69) 

(0.77) 

d Two-pulse mutated state with 26 roll pairs. 


hTwo-pulse mutated state with 24 roll pairs. 


c Two-pulse mutated state with 22 roll pairs. 



resulting values of a, b y v, and c g . There is no indication 
that a higher-order characteristic equation would lead to an 
improvement in these fits. 

In all cases the fitted values of v are positive, implying 
that the pulse (or pulse pair) travels more slowly than the 
group velocity c g for long-wavelength perturbations. More- 
over, in these cases the fitted values of a are positive, imply- 
ing that the basic wave train is unstable with respect to long- 
wavelength phase instabilities, as hypothesized. A significant 
source of error arises from fits that put too much weight on 
points that are contaminated by nonlinear effects; this prob- 
lem increases with decreasing interpulse distance and hence 
becomes particularly acute for multipulse states. In such 
states the double-pulse solution of the dynamical system (34) 
with the appropriate period may not approach the origin 
closely enough to sample the linear dispersion relation that 
describes its behavior in its vicinity. Consequently we expect 
that the values of the dispersion relation coefficients deduced 
from the corresponding fits will be less accurate than those 
obtained in the single-pulse case. This may account for some 
of the scatter in the coefficient values, particularly in a. 


TABLE VI. The values of the dispersion relation coefficients a , 
b, and v in Eq. (36) obtained from the fitted values of the eigen- 
values s. 


R r IR i “ ) 

a 

b 

V 

c * 

1.19 

0.42 

-1.90 

2.21 

8.09 

1.23 

0.91 

-2.21 

2.78 

8.20 

l.23 a 

1.79 

-2.78 

3.10 

8.54 

1.23 b 

0.22 

-2.35 

2.94 

8.31 

1.23 c 

0.61 

-2.00 

3.73 

9.39 

1.27 

0.21 

-1.60 

1.46 

5.95 

d Two-pulse mutated state with 26 roll pairs. 
hTwo-pulse mutated state with 24 roll pairs. 
c Two-pulse mutated state with 22 roll pairs. 




These results parallel a recent measurement of the phase dif- 
fusion coefficient for traveling wall modes in rotating con- 
vection [29]. See also [20]. 

The last column in Table VI lists the values of the group 
velocity deduced from the values of v in the previous col- 
umn with the help of the results listed in Table IV. The rough 
consistency of these values is reassuring, and supports our 
interpretation of the pulse formation process as the conse- 
quence of a long-wavelength phase instability and its de- 
scription using the third-order phase equation (29). Note that 
such consistency could not hold if we fitted each pulse sepa- 
rately. Such individual fits would yield a value of v but this 
value would reflect the speed of the pulse only in the absence 
of all others. In fact the pulse pairs travel together with a 
common (constant) speed, determined by the pulse-pulse in- 
teraction, and it is this speed that is determined by fitting the 
largest interpulse interval. This proviso applies to the single 
pulse as well; because of the use of periodic boundary con- 
ditions a single-pulse state is in fact a periodic train of 
equally spaced pulses. Such a pulse train has a different 
speed than an isolated pulse of identical shape. 

E. Pulse multiplicity and the ShiEnikov mechanism 

The inference that Eq. (29) possesses a homoclinic orbit 
of saddle-focus type for appropriate values of the coefficients 
of the nonlinear terms has a number of additional conse- 
quences. A general theory of the dynamics (in f) when such 
orbits are present is due to ShiEnikov [30]; its elaboration by 
Glendinning and Sparrow [31] is particularly useful. The 
theory defines a quantity S (the eigenvalue ratio) as 


and shows that a homoclinic orbit to a fixed point with S 
> 1 is generically simple and isolated. On the other hand, if 
£< 1, the orbit coexists with an uncountable number of non- 
periodic orbits as well as a countable number of periodic 
orbits of arbitrarily large periods. Although all these orbits 
correspond to different spatial states they do not resemble the 
observed pulselike states. However, Glendinning and Spar- 
row show that for parameters near those required for the 
formation of the primary (single-pulse) homoclinic orbit 
there are other parameter values (i.e., values of u) for which 
double-pulse homoclinic orbits are present. These ho- 
moclinic orbits (called secondary by Glendinning and Spar- 
row) start out as if to form a single pulse but miss the origin 
the first time around, and so spiral away from it again, form- 
ing a second pulse. Only after a second large-amplitude ex- 
cursion from the origin do they connect to the origin. Such 
solutions (see Fig. 11) are exactly of the form of our two- 
pulse states. Note that our approach predicts that a double 
pulse should move at constant speed given by the corre- 
sponding eigenvalue u, maintaining constant separation be- 
tween the two pulses. As suggested in Fig. 1 1 there is in fact 
a countably infinite number of such double-pulse (also triple- 
pulse etc.) states which accumulate (exponentially) on the 
parameter value, v* say, at which the primary one-pulse 
homoclinic orbit is present. These orbits differ in the number 
of turns each makes around the origin: the outer ones (v 
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FIG. 11. The modulation wavelength A vs v showing the cre- 
ation of single- and multiple-pulse states near the primary ho- 
moclinic connection at v* (corresponding to a single pulse on the 
real line), satisfying the Shil’nikov condition S<\. 


-y* large) perform few turns and hence correspond to 
double pulses with few oscillations in q separating them, 
while the inner ones perform more and more such oscilla- 
tions as their eigenvalue u approaches v* (see Fig. 11). The 
double-pulse homoclinic orbits accumulate on v* at a geo- 
metric rate given by, cf. [32], 


v n + \ 



2 TTfJi 

V 


+ 0(b) 





yja + 3/x 


+ 0(b) j, 
(42) 


where p* = 2/z(a + 4/x 2 ). In contrast the saddle-node bifur- 
cations (see Fig. 11) accumulate at the rate 


^ii = -exp ( -— + 0(b ) ) 

v„ \ V ! 

= ~ exp( - I == ^ + 0(fo)j. (43) 

\ \a + 3 /a I 

It should be noted that no similar behavior occurs in case (ii) 
in which only the primary one-pulse orbit is present at u*, 
cf. [33]; as v^v* the period of a periodic array of pulses 
increases monotonically with no horseshoe formation, 
much as in case (i) with S> 1 . 

Thus far we have discussed the consequences of having a 
homoclinic orbit to a saddle focus satisfying the Shil nikov 
condition 0< S< 1. This theory applies on the real line 
<£< + x. Our numerical simulations use a finite albeit large 
domain with periodic boundary conditions. Consequently, 
true homoclinic pulses cannot form, and our pulse states 
must be interpreted as periodic states with a long but finite 
period. In particular, the plethora of nonperiodic trains of 
pulses that coexist with the primary pulse at u* cannot be 
realized in our periodic geometry. Among the periodic states 
present at t? * only those with the correct period can be real- 
ized. We indicate this process by the dashed horizontal lines 


drawn in Fig. 1 1 which correspond to the period V and its 
integer fractions TIN, N = 1 , 2 ,... .As the wavelength 
curve A(v) approaches the primary homoclinicity at u* it is 
cut by the dashed lines. The intersection with the line TIN 
indicates the presence of a single-pulse state with period V /N 
at the corresponding v. This state can be replicated forming 
an A-pulse state that fits within the full period T . Such mul- 
tipulse states are characterized by equal separation between 
the constituent pulses. Figure 1 1 also shows one of infinitely 
many subsidiary homoclinicities corresponding to the forma- 
tion of double pulses. The speed of such pulses corresponds 
to the intersection of the corresponding A(v) loci with the 
line T. These double-pulse states are fundamentally different 
because while they still move rigidly with the speed given by 
this intersection their separation will be uneven. Moreover, 
as indicated in Fig. 1 1, the different double states formed this 
way will have different numbers of oscillations in their tails 
and will propagate with slightly different speeds than pairs of 
identical pulses. Figures 2(b) and 8(b) illustrate two distinct 
solutions, both with 26 roll pairs and two pulses with uneven 
pulse separation. In a similar way one can construct a variety 
of three-pulse states [see Fig. 8(a)]. However, recent theoret- 
ical work [34] indicates that such constructions are not of 
finite codimension; in particular, on the real line three-pulse 
states can be produced by arbitrarily small changes in the 
eigenvalue ratio S. Consequently any attempt to describe the 
observed three-pulse states by a truncated dynamical system 
of the form (34), i.e., a truncated phase equation, is doomed 
to fail. 

In the theory the case S= \ is a degenerate one because 
the linearization about the origin is divergence-free. This is 
the case for the Kuramoto-Sivashinsky equation [35]. In the 
present case when dispersion (b^O) is included we find that 

<44) 

Thus the ShiT nikov mechanism is generic in the present 
problem only in the presence of dispersion. 

F. Pulse solutions of the phase equation 

The conclusions we have drawn so far are predicated on 
the presence of a homoclinic orbit of saddle-focus type in 
Eq. (29) when the coefficients of the linear terms take the 
values deduced from the observed pulse shapes. But whether 
such an orbit is present is clearly a function of the nonlinear 
terms. It is a simple matter to show that for some choices of 
the coefficients the requisite homoclinic orbit cannot exist. 
For example, when u>0, b < 0 (and A=0) this is the case 
for all 0. Evidently a more detailed theory would start 
with known (calculated) values of the coefficients a,...,/t, 
and ask for what value (or values) of v do such orbits exist. 
Such a theory would not only predict the shape of the pulse 
but also its speed, rather than simply looking for consistency 
in the general picture. Since we are not in a position to com- 
plete this program we simply confirm that the phase equation 
(29) is capable of describing pulse states of the type de- 
scribed in Sec. III. To this end we have solved the nonlinear 
eigenvalue problem (33) on a finite domain of length I 
= 64 with a discretized one-dimensional mesh using an itera- 
tive Newton-Raphson-Kantorovich scheme [36] of high ac- 
curacy in both the Li norm of and the corresponding 
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R t = 1.19 R{."> 



FIG. 12. (a) Solution of the phase equation (33) with a 

= 0.4, b= — 1 .9, c= 1.0, e = — 800, g = 820, * = 700, and an im- 
posed phase jump A = 3.63X 2ir. The eigenvalue e = 6.20 implies a 
pulse speed u p “=1.88. (b) The instantaneous local wave number 
q(x) = 2ir/\(x) computed from Fig. 10(a) for comparison with (a). 

eigenvalue v . For the linear coefficients we used the values 
deduced from the pulse shape in Sec. IV D, a = 0.4, b~ 
~ 1.9, c= 1.0; in an effort to model the single pulse observed 
at R t = 1.19 /?j? f we used <p(0) = 0, <p(D=A, as the bound- 
ary conditions, with periodic and A determined from 

Fig. 10(a). The result shown in Fig. 12(a) is for e<0 and 
A' >0 for which homoclinic orbits to a nonzero background 
wave number are possible. The figure shows a single pulse in 
terms of the wave number modulation <pA£), riding on a 
sloping background corresponding to a nonuniform back- 
ground wave number. In Fig. 12(b) we show, for compari- 
son, Fig. 10(a) redrawn in terms of the local wave number. 
For the parameter values chosen, the solution of the phase 
equation yields a broader pulse than realized in the simula- 
tion. As discussed in Sec. IV D the wave number gradient 
was found to be an essential part of the fitting procedure [Eq. 
(40)]. The asymmetry of the pulse in Fig. 12(a) is similar to 
that of the observed pulse, while the corresponding eigen- 
value v = 6.20 yields a pulse speed v p that is about one-third 
of the observed one. 

The phase equation (29) exhibits a rich variety of pulse- 
like solutions which will be described elsewhere. The formu- 
lation of the pulse formation problem on a finite domain 
facilitates the location of such solutions, but requires the 
specification of the phase jump A across the domain. The 
finite size of the domain tends to suppress phase turbulence, 
and extends the parameter regime in which the phase equa- 
tion exhibits pulselike solutions. 

G. Stability considerations 

Thus far we have not addressed the stability properties of 
the various states identified above as solutions of the phase 
equation (29). We plan to give a detailed discussion of these 
and of the multiplicity of pulselike solutions to Eq. (29) else- 
where. However, the following remarks can be made. 


We have taken the point of view that the states found in 
the numerical simulations of Eqs. (l)-(5) are stable, al- 
though this is in fact something that is extremely hard to 
demonstrate, both mathematically and numerically (because 
of the very long horizontal diffusion times across the do- 
main). Any stability discussion divides naturally into two 
parts, the stability of a single-pulse state, and the stability of 
a periodic train of pulses. Assuming that only one pulse is 
present per period one first examines the stability of this state 
with respect to perturbations with the same wavelength; sub- 
sequently one can examine the stability with respect to per- 
turbations with periods that are increasing multiples of the 
basic period. Unfortunately the stability properties of ho- 
moclinic pulses on the whole real line are not easily estab- 
lished, even in dissipative systems [37], In fact we expect 
that most such states will prove to be unstable, possibly in- 
cluding the one shown in Fig. 12(a). Usually stable pulses on 
the real line require that the equilibrium state that is con- 
nected by the pulse is stable; otherwise the pulse inherits the 
instability of the equilibrium. Thus stable pulses are usually 
subcritical. In our case this argument suggests that stable 
pulses should only be observed when the phase diffusion 
coefficient -a is positive. While this is so for the Eckhaus 
instability, in dispersive systems, such as the one considered 
here, stable pulses can be present even when the basic equi- 
librium is unstable, provided that the instability is only con- 
vective. In fact the theory described above applies equally 
well for a negative as for a positive. The only difference is 
that when a<0 the quantity fi in Eq. (37) must be suffi- 
ciently large. It follows that v must also be sufficiently large 
and hence that the pulse speed v p will differ substantially 
from the group velocity c g . In fact this is the case for the 
pulses we have found (see Table IV). However, we have 
been unable to find any signature of the required hysteresis 
in the pulse formation process and consequently favor the 
choice a>0. Assuming that a single pulse per wavelength is 
stable the stability of a whole train of such pulses may then 
be described as in Refs. [38,39], One finds that instability 
that leads to pulse bunching (i.e., a nonuniform train of 
pulses) depends on the dispersion relation giving the pulse 
speed v p as a function of the basic period or wavelength. 
Unfortunately since all our calculations were done with the 
same aspect ratio T = 64 we have no information about this 
dispersion relation and hence about the stability properties of 
whole trains of pulses. 

Finally, it is of interest to note that the phase equation 
examined by Janiaud et al. [20] also exhibits the formation 
of pulses, although in a rather limited parameter regime lo- 
cated astride the Benjamin-Feir instability boundary and very 
close to it. In the equation studied here the pulse formation 
process appears to be more robust, while its special form 
facilitates analytical approach. 

V. DISCUSSION AND CONCLUSIONS 

In this paper we have initiated a study of the mechanisms 
that lead to the formation of confined states in two- 
dimensional large-aspect-ratio thermosolutal convection. Nu- 
merical simulation of the Boussinesq equations describing 
the system in a two-dimensional periodic domain with real- 
istic boundary conditions revealed the existence of a plethora 
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of confined states, both single pulse and multipulse, all ot 
which were found to travel rigidly on top of a uniform train 
of fully nonlinear traveling waves. All of our pulses traveled 
in the same direction as the underlying wave train but did so 
more slowly, sometimes substantially so, and appeared to be 
created spontaneously as a result of a phase instability. In 
particular, no finite-amplitude perturbations were required to 
form these pulses, and no evidence of hysteresis in their 
formation was uncovered. The pulses lead to the appearance 
of spatially nonuniform wave trains in which patterns with 
differing phase velocities, amplitudes, and wavelengths co- 
exist stably. They are characterized by a larger degree of 
coherence between the vertical velocity and the thermal and 
solutal fields and as a result transport more heat and solute 
than the background state. Consequently, the process of 
pulse formation can be interpreted as an attempt by the sys- 
tem to maximize heat and solute transport within the con- 
straint imposed by phase conservation, i.e., a fixed number of 
rolls. 

We have focused on the Rayleigh number regime 

1.15/?j^ > </?r < 1.27/?V ?) within which P ulse formation is ini " 

tiated but the phase of the pattern remains well defined. This 
regime allows one to employ the techniques of phase dynam- 
ics, one of the few analytical means available to us for un- 
derstanding confined states that develop on top of a highly 
nonlinear state. We have used this approach to gain substan- 
tial insight into the properties of these pulses and their mul- 
tiplicity. These pulses appear to be well described by a third- 
order phase equation. This equation is a particularly simple 
case of the phase equation describing the evolution of long- 
wavelength Benjamin-Feir instability of traveling waves 
[20]. To our knowledge this important special case has not 
hitherto been studied. We have found, using our phase equa- 
tion, that the pulse shape is given by a third-order nonlinear 
eigenvalue problem, with the eigenvalue determining the 
pulse speed. The theory makes predictions about the pulse 
shape that are consistent with the results of our simulations. 
These include the fact that pulses traveling more slowly than 
the group velocity should have oscillatory tails, and that the 
local wavelength and phase velocity should be linearly re- 
lated. In addition we verified a general expression relating 
the pulse speed to the wavelengths and phase velocities 
within the pulses and without. We interpreted the oscillatory 
tail of the pulses in terms of a homoclinic orbit to a saddle 
focus satisfying, under mild assumptions on the phase diffu- 
sion and dispersion coefficients, the ShiFnikov condition 0 
<S<\. Using standard theory we established the presence, 
on the real line, of a great number of coexisting pulse and 
nonpulse states, and applied this result to establish the exis- 
tence on a periodic domain of not only single pulses but also 
of two different types of double (and triple etc.) pulses, all of 
which move rigidly at constant speed given by the eigen- 
value of an appropriate nonlinear eigenvalue problem. These 
included periodic two-pulse states in which the pulses are 
identical and their separations equal, and others in which 
they are not. We exhibited examples of such pulses in the 
numerical simulations, obtaining them by a process we 
called 4 ‘mutation” in which selected parts of the pulse were 
culled and replicated to create an initial condition of the de- 
sired kind. Time integration forward in time was then used to 
check their stability. 


Despite the success of our approach we caution that the 
phase approach has an important shortcoming. Formally it 
applies only to states with weak phase gradients. In our 
pulses wave number changes of order one were common. 
Thus our states are outside the regime in which the equation 
should hold. It appears, however, that phase conservation is a 
sufficiently strong constraint that the phase equation provides 
a reliable model even outside its formal range of validity. 

Another approach to the present problem would be to 
consider a pulse as the bound state of two fronts, with the 
pulse width determined by the minimum of the potential be- 
tween them (if one exists). The small oscillations observed in 
the widths of the pulses during their propagation (see Sec. 
IV) could then be interpreted as oscillations about the mini- 
mum of this potential. Although such oscillations should de- 
cay with time, we have found them to be long lived. Al- 
though the most likely explanation of such oscillations 
(based on their frequency) is that they are an effect of the 
finite size of the rolls we cannot exclude the possibility that 
they are the consequence of self-excited oscillations of such 
bound states (cf. [40]). This approach, which offers scope for 
a more realistic model of the confined states described here, 
is planned for a future paper. 

The present paper represents the first observation of con- 
fined states in thermosolutal convection. These confined 
states are distinct from the ones found in binary fluid con- 
vection. The latter occur very close to onset and consist of 
localized packets of traveling waves on an otherwise undis- 
turbed background. They also differ physically. In binary 
fluid mixtures there is no net concentration transport, only 
concentration migration in response to the applied tempera- 
ture gradient. Thus the concentration is ultimately pumped 
horizontally, out of the packet, where it creates regions of 
negative buoyancy. As a result the packets travel much more 
slowly than might otherwise be expected. However, because 
these states are found near onset they are still accessible to 
description by weakly nonlinear theory [5]. In contrast, in the 
thermosolutal system solute is transported in the vertical. As 
a result there is very little tendency towards the kind of con- 
centration braking that is so important in binary fluids. Re- 
cent developments in experimental technique [8] suggest that 
the pulses described here should be amenable to experimen- 
tal study. The physical considerations just alluded to also 
suggest that any small-amplitude subcritical localized states 
present near onset of oscillatory thermosolutal convection 
should travel faster and have different stability properties 
than the corresponding ones in binary fluid convection. How- 
ever, for our parameter values there is an additional differ- 
ence as well. The subcritical traveling waves present at these 
parameter values have a second unstable eigenvalue, they are 
unstable with respect to standing waves, which also bifurcate 
subcritically. Consequently we do not expect localized pack- 
ets of traveling waves near onset in our system, although the 
possibility of small-amplitude confined states of standing 
waves cannot be excluded. 

In their study of localized traveling waves in binary fluids 
near onset Barten et al. [7] showed that the localized states 
drifted with a small group velocity in the propagation direc- 
tion of the phase of its TW components; in the comoving 
frame these states were found to be time periodic. Different 
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localized states differed only in the width of the central part 
of the pulse while the leading and trailing ends (although 
different) were found to be the same in all the simulations 
performed. The velocity, temperature, and concentration pro- 
files could not be described by one common amplitude since 
their widths and shapes were different. The velocity and tem- 
perature amplitudes in the center of the localized states were 
only slightly smaller than in the extended states at the same 
values of the parameters. However, the local wavelengths 
and phase speeds of the rolls increased monotonically as they 
moved from the trailing toward the leading end. The fre- 
quency of the traveling waves in such a pulse was found to 
be about half the Hopf frequency at the onset, a much larger 
value than in an extended state, indicating a larger concen- 
tration contrast between adjacent rolls and a greater mean 
lateral concentration flux for the localized states than for the 
extended ones. The wavelength in the pulse was only about 
10% smaller than in an extended state of similar frequency. 

Barten et al. found that the mean buoyancy force in such 
convection is dominated by the concentration contribution 
rather than the thermal one. A phase shift between the con- 
centration wave and the velocity wave which takes place in 
the central part of the localized state drove a mean concen- 
tration circulation extending over the whole localized state. 
This in turn induced a large-scale concentration redistribu- 
tion. Ahead of the leading front this concentration redistribu- 
tion produces a barrier that weakens the mean buoyancy 
force, thereby impeding a rapid invasion of the conductive 
region by convection. On the other hand, the different widths 
of the concentration and temperature pulses appear to en- 
hance the buoyancy under both fronts and thereby stabilize 
the localized state against invasion by the conduction state. 
The mean concentration flux decreases with decreasing fre- 
quency (increasing R r ) and with it the concentration redis- 
tribution. This weakens the barrier and increases the forward 
drift motion of the localized state. Larger frequencies 
strengthen the barrier and the forward drift velocity of the 
pulse into the quiescent regime is reduced. These results are 
in good agreement with the experimental results of Kolodner 
[41] and with the analytical theory of Riecke [5]. We have 
identified similar large-scale communication in thermo- 
solutal convection using moving averages. In the reference 
frame moving with a pulse the small-scale dynamics are fil- 


tered out, revealing the dynamics on the scale spanning the 
full horizontal width of the domain and involving large-scale 
vortices that drift with the pulses. The pulses owe their ex- 
istence to an interaction between the order one convective 
scale and this large scale. 

We should remark that traveling wave thermosolutal con- 
vection has not hitherto been simulated with no-slip bound- 
ary conditions, even in domains with basic period of order 
the roll wavelength. Here, however, our results do not differ 
markedly from those obtained by Barten et al. for binary 
mixtures [16]. In particular, the velocity and temperature 
fields are observed to have similar structure and strength to 
those in comparable stationary states in pure-fluid convec- 
tion, with the temperature wave phase shifted slightly rela- 
tive to the velocity wave, cf. [1]. In contrast to the velocity 
and temperature components, which at our Rayleigh num- 
bers were nearly harmonic, the solutal component was found 
to be highly nonharmonic, with an almost trapezoidal profile, 
indicative of nearly homogeneous solute within each roll, 
separated by a thin boundary layer that snakes alternately 
between the top and bottom of the layer. This is a conse- 
quence of the “mirror-glide” symmetry which in turn is a 
consequence of imposing identical boundary conditions at 
the top and bottom of the layer. The lateral mean flow in a 
uniform wave train was found to be several orders of mag- 
nitude smaller than the phase velocity of the wave. A small 
but significant mean lateral convective heat flux, ranging be- 
tween 0.006 N T and 0.014 N r , is generated by the phase 
shift between the vertical velocity of the fluid and the tem- 
perature wave, with a similar solute flux, between 0.01 5 N s 
and 0.038A 5 , generated by the phase difference between the 
convection and the solutal wave. Such mean flows and fluxes 
were also found in simulations of binary fluid mixtures [16]. 
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